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ABSTRACT 

We propose to use high-redshift long 7-ray bursts (GRBs) as cosmological tools to constrain 
the amount of primordial non-Gaussianity in the density field. By using numerical, N-body, 
■ hydrodynamic, chemistry simulations of different cosmological volumes with various Gaus- 

fH I sian and non-Gaussian models, we self-consistently relate the cosmic star formation rate den- 

O ■ sity to the corresponding GRB rate. Assuming that GRBs are fair tracers of cosmic star for- 

mation, we find that positive local non-Gaussianities, described in terms of the non-linear 
parameter, /nl, might boost significantly the GRB rate at high redshift, z ^> 6. Deviations 
with respect to the Gaussian case account for a few orders of magnitude if /nl~ 1000, one 
order of magnitude for /nl^ 100, and a factor of ~ 2 for /nl^ 50. These differences are 
found only at large redshift, while at later times the rates tend to converge. Furthermore, a 
comparison between our predictions and the observed GRB data at z > 6 allows to exclude 
large negative /nl, consistently with previous works. Future detections of any long GRB at 
extremely high redshift (z ~ 15 — 20) could favor non-Gaussian scenarios with positive /nl- 
More stringent constraints require much larger high-z GRB complete samples, currently not 
available in literature. By distinguishing the contributions to the GRB rate from the metal- 
^1 ■ poor population III regime, and the metal-enriched population II-I regime, we conclude that 

the latter is a more solid tracer of the underlying matter distribution, while the former is 
strongly dominated by feedback mechanisms from the first, massive, short-lived stars, rather 
than by possible non-Gaussian fluctuations. This holds quite independently of the assumed 
population III initial mass function. 
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1 INTRODUCTION for which the cosmological constant A represents the simplest ex- 

planation. Thanks to the evidences coming from different observa- 
The present standard cosmological model assumes that a pnmor- , , . . , . , .. ,, ji 

r a r tional datasets (mainly cosmic microwave background, galaxy sur- 

dial inflationary phase (Starobinsky 1980; Guth 1981; Linde 1990) , , f TT . 

\ „ . veys and supernovae), the general properties 01 our Universe have 

ends with the creation of density fluctuations, that then grow during , , . .. , . , , , 

become clearer and its parameters known with much better accu- 
cosmological times (e.g. Gunn & Gott 1972; Weinberg 1972; Press ™ . , . ■ , v 

racy. The estimated contributions to the cosmic density are (Ko- 

fi Schechter 1974; White & Rees 1978; Peebles, P. J. E. 1993; 

Sheth & Tormen 1999; Peacock, J. A. 1999; Hogg 1999; Coles, 



matsu et al. 2011): ft ,m = 0.272, Q. , A = 0.728, Q. , b = 0.044, 

for matter, cosmological constant, and baryons, respectively; the 
P. & Lucchin, F. 2002; Peebles & Ratra 2003) to give birth to the . ' . , ° t ' . . { ' F V 

cosmic equation 01 state parameter is consistent with w — — 1, the 



presently observed large scale structure of the Universe (Barkana 
& Loeb 2001; Ciardi & Ferrara 2005; Bromm & Yoshida 2011) 



theoretical expectation of the cosmological constant; the primor- 
dial power spectrum has spectral index n = 0.96, and a normaliza- 
Stars, galaxies, and clusters of galaxies form by gravitational col- .. ,. . /, . , 

, a of n corresponding to a mass variance within a 8 Mpc//i- sphere of 



lapse in an expanding flat space, composed by ~ 30% of matter 
and ~ 70% of an unknown constituent referred to as dark energy, 
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a 8 = 0.8. 

Even if the above picture is quite satisfying, the specific mecha- 
nism driving the inflation is however not completely understood. 
This fact justifies the existence in the literature of a plethora of 
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possible inflationary models, each of them with specific predic- 
tions for various observables. In particular, the study of the statis- 
tical distribution of the primordial fluctuations is considered one 
of the best way to discriminate between them. In fact, alterna- 
tives to the standard single-field slow-roll model, which predicts 
a nearly Gaussian distribution, can produce a significant amount of 
non-Gaussianity (Bartolo et al. 2004; Chen 2010). The most recent 
analyses of the observational data show some evidence for possible 
departures from Gaussianity, even if with a low level of signifi- 
cance (Peebles 1983; Desjacques & Seljak 2010; Lo Verde & Smith 
2011; D'Amico et al. 2011; Komatsu et al. 2011; as also collected 
in the summary Table 2 by Maio & Iannuzzi 201 1). More precisely, 
slightly positively skewed models are favored. 
From a theoretical point of view, the presence of some amount 
of primordial non-Gaussianity has two main effects, which are 
then used as efficient constraining tools: it introduces a scale- 
dependence in the bias factor (e.g. Grinstein & Wise 1986; Mc- 
Donald 2008; Desjacques et al. 2009; Fedeli et al. 2011; Norena 
et al. 2012; Wagner & Verde 2012), and it modifies the abundance 
and the formation history of rare events (i.e. very low- and high- 
sigma fluctuations; e.g. Koyama et al. 1999; Zaldarriaga 2000; 
Grossi et al. 2007, 2008; Wagner et al. 2010; LoVerde & Smith 

2011) . High redshifts represent an interesting regime to potentially 
investigate these effects. Indeed, very early structures and primor- 
dial mini-haloes hosting the first bursts of star formation are ex- 
pected to be somehow affected by the presence of primordial non- 
Gaussianities (as discussed in Maio 201 1). 

In more detail, due to the sensitivity of the gas cooling capa- 
bilities to the underlying matter density field, numerical hydro- 
dynamical simulations have shown that the initially skewed non- 
Gaussian features could be reflected by the probability distribu- 
tion function of the high-z cosmic medium (Viel et al. 2009), by 
a change in the molecular gas evolution and formation epoch of 
first stars and galaxies (Maio & Iannuzzi 201 1; Maio 201 1), and by 
the consequent metal pollution in the Universe (Maio & Khochfar 

2012) . Furthermore, simple semi-analytical arguments have sug- 
gested non-Gaussianity effects on the birth of primordial black 
holes (e.g. Bullock & Primack 1997; Green & Liddle 1997; Ivanov 
1998; Avelino 2005; Hidalgo 2007; Kohri et al. 2008; Bugaev & 
Klimai 2012; Byrnes et al. 2012), cosmic reionization (Crociani 
et al. 2009), and hydrogen 21-cm signal (Cooray 2006; Cooray 
et al. 2008; Pillepich et al. 2007; Joudaki et al. 2011; Chongchit- 
nan& Silk 2012). 

In this respect, a key tool for studies of high-redshift environments 
might be represented by 7-ray bursts (GRBs), powerful explosions 
emitting 7 rays in the ~ [1 keV, 10 MeV] energy band, mostly 
distributed around ~ 0.1 — 1 MeV, as detected by the latest Fermi- 
GBM instrument operating in the [8 keV, 40 MeV] range (Bissaldi 
etal. 2011). 

These bursts have: isotropic equivalent peak luminosities as high 
as ~ 10 54 erg s _1 (the record holder being GRB 080607, Per- 
ley et al. 2011); an isotropic angular distribution (Fishman et al. 
1994; Paciesas et al. 1999, 2012); and a bimodal duration distri- 
bution (Kouveliotou et al. 1993), with most of them lasting for a 
period longer than 2 seconds (long GRBs), and some of them, de- 
tected mostly at low redshift, for a period shorter than 2 seconds 
(short GRBs). 

In the following we will only consider long GRBs (LGRBs), which 
are supposedly related to the death of massive stars (see extensive 
reviews by e.g. Piran 2004; Meszaros 2006) and, therefore, they are 
indicators of the local star formation episodes (e.g. Jakobsson et al. 
2005; Nuza et al. 2007; Lapi et al. 2008; Yiiksel et al. 2008; Kistler 



et al. 2009; Butler et al. 2010; Campisi et al. 201 1; Mannucci et al. 
2011; Ishida et al. 201 1; Elliott et al. 2012; Robertson & Ellis 2012; 
Michalowski et al. 2012). 

Typical 7-ray bursts have long-lasting afterglows at lower frequen- 
cies, from the X-rays to the radio band due to scattering with 
the surrounding ambient medium (Paczynski 1991; Dermer 1992), 
and are theoretically explained by the "collapsar model" (Woosley 
1993; Wang & Wheeler 1998; Meszaros et al. 1999; Woosley & 
Heger 2012): a massive black-hole stellar remnant - probably a 
Wolf-Rayet star (but see Baron 1992; Yoon et al. 2010) - accret- 
ing stellar mass from a disk (Popham et al. 1999; Fryer et al. 1999; 
Narayan et al. 2001; Yoon & Langer 2005; De Colle et al. 2012) 
at a rate of ~ 0.01 — 10 M /s, and accompanied by a collimated- 
jet emission with a few degree opening angle (e.g. Waxman 1997; 
Rhoads 1997; Sari et al. 1998; Wang & Wheeler 1998; Schmidt 
1999, 2001). 

Due to additional factors, like asymmetric explosions or stellar ro- 
tation (e.g. Sollerman et al. 2005; Thone et al. 2008; Ostlin et al. 
2008), only a small fraction of SNe, ~ 10~ 2 - 10~ 3 (Fruchteret al. 
2006; Yoon et al. 2006; Bissaldi et al. 2007; Soderberg et al. 2010; 
Grieco et al. 2012), can result into a LGRB. However, also taking 
into account such effects, there is still significant lack of knowledge 
of some important details, like the minimum mass for GRB black- 
hole progenitors, that is highly debated and expected to lie between 
typical SN limits of ~ 25 — 40 Mq (see also recently proposed 
upper values of even ~ 60 M Q in Georgy et al. 2012). 
The uniquely bright luminosities of GRBs facilitate their detection 
up to very high redshift, as shown by the three bursts spectroscop- 
ically confirmed at z > 6, i.e. GRB 050904 at z = 6.3 (Kawai 
et al. 2006), GRB 080913 at z = 6.7 (Greiner et al. 2009), and 
GRB 090423 at z = 8.2 (Salvaterra et al. 2009; Tanvir et al. 2009), 
and by the case of GRB 090429B, having a photometric redshift of 
z ~ 9.4 (Cucchiara et al. 201 1). 

High-redshift GRBs are a powerful and, in some cases, a unique 
tool to study the Universe at the early stages of structure formation 
and can provide fundamental information about the environment of 
their own hosting galaxies like: 

• metallicity and dust content (Savaglio et al. 2005; Savaglio 
2006; Nuza et al. 2007; Fynbo et al. 2008; Savaglio et al. 2009; 
Mannucci et al. 2011; Campisi et al. 2011; Niino et al. 2011); 

• neutral-hydrogen fraction (Nagamine et al. 2008; McQuinn 
et al. 2008; McQuinn & et al. 2009; Gallerani et al. 2008; Mirabel 
et al. 2011; Robertson & Ellis 2012); 

• local inter-galactic radiation field (Inoue et al. 2010); 

• early cosmic magnetic fields (Takahashi et al. 201 1); 

• stellar populations (Bromm & Loeb 2006; Campisi et al. 20 1 1 ; 
de Souza et al. 2011; Salvaterra et al. 2012); 

In the present work, we argue that GRBs can be additionally used 
as cosmological probe of the amount of non-Gaussianity present in 
the primordial density field. In fact, they are sensitive to the under- 
lying cosmological model through the first episodes of the cosmic 
star formation history. 

We will show how GRBs trace the matter distribution at high red- 
shift by performing a detailed analysis of the GRB rate in different 
non-Gaussian scenarios, with the help of N-body, hydrodynamic, 
chemistry simulations of early structure formation (Maio & Ian- 
nuzzi 201 1). In the simulated volumes, star formation is addressed 
on the basis of the local thermodynamical properties of the col- 
lapsing gas, by consistently following its density, temperature and 
chemical composition, and by taking into account stellar evolution 
and feedback effects. The resulting star formation rate (SFR) and 
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the adopted initial mass function (IMF) for the stellar populations 
tracked during the runs are used to determine the expected GRB 
formation rate density in the various cases, and hence the inte- 
grated GRB rate (R), for both metal-poor population III (hereafter, 
popIII) regime and metal-enriched population II-I (hereafter, popll- 
I) regime. 

The paper is structured as follows. In Sect. 2 we describe the 
numerical simulations used in our study; in Sect. 3 we present 
the adopted model for GRB evolution (Sect. 3.1), its valida- 
tion (Sect. 3.2), and the consequencies for non-Gaussian models 
(Sect. 3.3); finally, in Sect. 4 we discuss and summarize our find- 
ings and conclude. In the following, when mentioning GRBs we 
will refer to LGRBs. 



2 NUMERICAL SIMULATIONS 

In the present paper, we will consider a set of numerical N- 
body, hydrodynamical, chemistry simulations with two different 
box sizes starting from initial conditions with a different level of 
primordial non-Gaussianity. A more detailed description of the 
simulations can be found in Maio & Iannuzzi (2011). Local non- 
Gaussianities were included in the initial conditions by adding 
second-order perturbations to the Bardeen gauge-invariant poten- 
tial (e.g. Salopek & Bond 1990): 

$ = $l + /nl [$l- < *l >] , (1) 

where 3>l is the linear Gaussian part, and /nl the dimensionless 
coupling constant controlling the magnitude of the deviations from 
Gaussianity in the large-scale-structure formalism. 
The simulations were performed by using a modified version of 
the parallel tree-PM/SPH Gadget code (Springel 2005), which in- 
cludes gravity and hydrodynamics, with radiative gas cooling both 
from molecular and atomic (resonant and fine-structure) transitions 
(Maio et al. 2007), a multi-phase model for star formation (Springel 
& Hernquist 2003), UV background radiation (Haardt & Madau 
1996), wind feedback (Springel & Hernquist 2003; Aguirre et al. 
2001), chemical network for e~, H, H+, FT, He, He+, He ++ , 
H 2 , Hj, D, D+, HD, HeH + (e.g. Yoshida et al. 2003; Maio et al. 
2006, 2007, 2009; Maio 2009; Maio et al. 2010, and references 
therein), and metal (C, O, Mg, S, Si, Fe) pollution from popIII 
and/or popII-I stellar generations, ruled by a critical metallicity 
threshold of Z crit = 10~ 4 Z© (see discussion in Tornatore et al. 
2007; Maio et al. 2010, 2011). The cosmological parameters are 
fixed by assuming a flat concordance ACDM model with matter 
density parameter fio, m = 0.3, cosmological-constant density pa- 
rameter S1o,a = 0.7, and baryon density parameter fi ,b = 0.04; 
the present Hubble parameter is fixed to Ho = 100 h km/s/Mpc, 
with h = 0.7. Finally, the matter power spectrum has a spectral 
index n = 1 and is normalized assuming that the mass variance 
within 8 Mpc/ft-radius sphere is erg = 0.9. 

A Salpeter IMF with mass range [0.1, 100] M© was adopted for 
the popII-I regime, while a top-heavy IMF with short-lived stars in 
the mass range [100, 500] M© was assumed for the popIII regime 
(see literature for further studies on the expected range of mas- 
sive popIII stars: Abel et al. 2002; Yoshida et al. 2003; Inayoshi 
& Omukai 2012; or low-mass popIII stars: Yoshida 2006; Yoshida 
et al. 2007; Campbell & Lattanzio 2008; Suda & Fujimoto 2010; 
and the impacts of the different assumptions: Maio et al. 2010). 
Massive stars die as SN or as pair-instability SN (PISN) in the 
range [8, 40] M© and [140, 260] M©, respectively, polluting the 



surrounding medium and enhancing the transition from a metal- 
poor to a metal-rich regime (e.g. Tornatore et al. 2007; Maio et al. 
2010, 2011). Black-hole remnants form from stellar masses in the 
ranges [40, 100] M© (popII-I progenitors), [100, 140] M© (popIII 
progenitors), and [260, 500] M© (popIII progenitors). 
To follow with sufficient accuracy all the relevant scales at the dif- 
ferent cosmological epochs, we consider two sets of simulations. 
The first one assumes small boxes with side of 0.5 Mpc/h, and 
allows us to resolve the gas behaviour down to ~ pc scales at 
z ~ 9 — 30 (Maio & Khochfar 2012), with gas and dark-matter 
mass resolutions of 42.35 M© /h and 275.28 M© / h, respectively, 
and comoving softening of 4 pc/h. 

The box size of the second set is much larger, 100 Mpc/h, 
so that we can resolve galactic ~ kpc scales at lower redshift 
(Maio et al. 2011), since gas and dark-matter mass resolutions are 
3.39 x 10 8 M©//i and 2.20 x 10 9 M©/ft, respectively, and the 
comoving softening is 7.8 kpc/h. 

For both sets of simulations, different levels of primordial non- 
Gaussianity have been considered, namely /nl = 0, 10, 50, 100, 
and 1000. We highlight that current data seem to suggest positive 
/nl values, between and 100 (e.g. Komatsu et al. 2011), but in 
the present work we will consider the / NL =1000 case as well, as an 
extreme example. For further details we refer to Maio & Iannuzzi 
(2011). 

The star formation rate for both stellar population regimes extracted 
from these simulations will represent the fundamental input for our 
estimates of the GRB rates, as described in the following sections. 
For the sake of clarity, in Fig. 1 we re-propose the redshift evolu- 
tion of the star formation rate densities derived from our ten runs 
and widely discussed in Maio & Iannuzzi (201 1). These curves are 
the starting point of our following analyses. 



3 CALCULATION OF THE GRB RATES 

In the following section we will present the results for the GRB 
rates expected from our simulations. Our starting point is the co- 
moving SFR density, p*, tracked by the different runs as a function 
of z (Maio et al. 2010; Maio & Iannuzzi 2011; Maio 2011), from 
which we compute the evolution of the GRB formation rate den- 
sity, wgrb, and hence the corresponding integrated GRB rate, R. 
We will proceed as follows: in the first place, we will present a phe- 
nomenological model describing the redshift evolution of GRBs as 
observed by Swift (Sect. 3.1); then we will validate it against obser- 
vational data at z ^ 6 (Sect. 3.2), i.e. the epoch when the effects of 
primordial non-Gaussianities are expected to play a major role; and 
eventually we will apply it to an ideal instrument that is assumed to 
detect all the GRBs produced in the different cosmological scenar- 
ios (Sect. 3.3). 

3.1 Model description 

The basic features of the model are presented in Sect. 3.1.1, fol- 
lowed by the derivation of the best-fitting values for the model 
free parameters in Sect. 3.1.2. We stress that the parameters of the 
model are dependent on the whole cosmic star formation history, 
and, therefore, they do depend on the /nl values, too. 

3.1.1 Formalism 

The expected redshift distribution of "observed" GRBs can be com- 
puted once the GRB luminosity function (LF) and the GRB forma- 
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Figure 1. Top panels: total and popIII star formation rate density evolution as a function of redshift, z, for the small (left) and the large (right) boxes. Different 
lines refer to /jml = (solid black lines), 10 (dotted cyan lines), 50 (short-dashed green lines), 100 (long-dashed red lines), and 1000 (dot-dashed blue lines). 
Bottom panels: ratio between the popIII and the total star formation rate. The plots are taken from Maio & Iannuzzi (201 1). 



tion history have been specified (e.g. Porciani & Madau 2001; Fir- 
mani et al. 2004; Natarajan et al. 2005; Guetta et al. 2005; Daigne 
et al. 2006; Salvaterra & Chincarini 2007; Salvaterra et al. 2009; 
Dai 2009; Wanderman & Piran 2010; Cao et al. 2011; Salvaterra 
et al. 2012). 

We briefly recap here the adopted formalism and refer the interested 
reader to Salvaterra & Chincarini (2007); Salvaterra et al. (2009, 
2012) for more details. 

The observed peak photon flux, P, emitted by an isotropically ra- 
diating source at redshift z and corresponding luminosity distance 
cLl{z), as detected in the energy band E m i n < E < -E max , is 



P = 



(1 + *) 



S(E)dE, 



(2) 



where S(E) is the differential rest-frame photon luminosity of the 
source. To describe the typical burst spectrum we adopt a "Band" 
function with low- and high-energy spectral indices equal to —1 
and —2.25, respectively (see also Band et al. 1993; Preece et al. 
2000; Kaneko et al. 2006). 

The spectrum normalization is obtained by imposing that the 
isotropic-equivalent peak luminosity is 



L = 



/■lOMeV 

/ E S{E) dE. 



(3) 



To estimate the peak energy of the spectrum, E p , for a given L, we 
correlate E p and L as done in Yonetoku et al. (2004); Nava et al. 
(2012); Ghirlanda et al. (2012). 

Given a normalized GRB LF, ip(L), the observed number rate of 
bursts (in yr^ 1 ) at redshift z with peak photon flux, P, between Pi 
and P2 is 



N(z) 



dNp 1<P< p 2 (z) 
dt 



dz 



, dV(z') h GRB {z') 
dz' 



(l + O 

-Lp 2 (z 1 ) 

I ip{L')dL', 

lL Pl (z>) 



(4) 



where the factor (1+z') 1 accounts for cosmological time dilation, 



dV(z) 
dz 



d£l d 2 c (z) j^- 



(5) 



is the comoving volume element, dSl is the solid angle d c (z) is 
the comoving distance, H{z) is the expansion parameter (for more 
explicit details see e.g. Weinberg 1972; Hogg 1999), c is the speed 
of light, and riGRB(z) is the actual comoving GRB formation rate 
density as a function of redshift. 

Here, we assume that GRBs are good tracers of star formation, and 
thus that the GRB formation rate density is directly proportional to 
the SFR density (see further discussion in Sect. 4), i.e. 



«grb(z) = kp i ,(z), 



(6) 



where the normalization constant, k (whose dimensions are the in- 
verse of a mass), incorporates further not-well-known effects, like 
GRB beaming (Frail et al. 2001; Panaitescu & Kumar 2001; Rossi 
et al. 2002; Ghirlanda et al. 2007), efficiencies (Fruchter et al. 2006; 
Yoon et al. 2006; Bissaldi et al. 2007; Soderberg et al. 2010; Grieco 
et al. 2012), and black-hole production probability (depending on 
the adopted IMF). 

We will adopt (see next section for more details) a normalized GRB 
LF described by a single power-law with slope £ and decreasing ex- 
ponentially below a characteristic luminosity, L«, 



L 



(7) 
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Figure 2. Redshift evolution of the observed number rate of GRBs, N, 
at the sensitivity of Swift instrument, corresponding to a peak flux of 0.4 
ph s — 1 cm -2 in the [15, 150] keV band, and for the Swift field of view 
Af2 s = 1.4 sr. Model results for fj^L = and fNL=1000 are shown as 
dark orange and light blue shaded regions, respectively, taking into account 
the eiTors on the evolution parameter. The trends (without errors) for the 
models with /nl = 10, 50, 100 are shown with long-dashed, short-dashed, 
and dotted lines, respectively. The arrow refers to the lower limit on the 
rate of GRBs at z > 6, imposed by the four confirmed detections at these 
redshifts (see text). 



Then, we consider the possibility that the GRB LF evolves by set- 
ting L*(z) = I/*,o(l + z) 5 , where L*,o is the characteristic lumi- 
nosity at z — 0, and 8 is the evolution parameter. 
For simplicity, the normalization of tp(L) is included in k, and it is 
fixed when the GRB number rate in equation (4) is normalized to 
the rate observed at z = 0. 

From the previous relations we can finally compute the GRB rate 
(in units of yr -1 sr -1 ), 7?, as: 



R(z) 



dN(z) 

do. 



(8) 



i.e. by taking the derivative with respect to the solid angle of the 
GRB number rate in equation (4). 

3.1.2 Parameter estimation 

The values of the free parameters of the model (i.e. L+,o, £, k and 
8) are optimized separately for all the models, by using the SFRs 
obtained from the different cases in the 100 Mpc//j-size boxes. 
We proceed by minimizing the C-stat function (Cash 1979), jointly 
fitting the observed differential number counts in the [50, 300] keV 
band of BATSE (Stern et al. 2001) and the observed redshift distri- 
bution of bursts in a redshift complete subsample between z = 0.13 
and z = 5.47 of Swift bursts with photon fluxes in excess of 
2.6 ph s -1 cm -2 in the Swift [15, 150] keV band (for more de- 
tails see Salvaterra et al. 2012). 

While the redshift complete Swift subsample provides a powerful 
test for the existence and the redshift evolution of the long GRB 
population, the fit to the BATSE number counts allows to obtain the 



normalization k and to better constrain the GRB LF free parame- 
ters. It is worth noting that the same best-fit parameters provide a 
good fit also to the Swift differential peak-flux number counts once 
the energy band ([15, 150] keV), the field of view (Af2 s = 1.4 sr), 
and the observing lifetime of Swift are considered (Salvaterra & 
Chincarini 2007). 

The best-fit values together with their 1-cr confidence levels are pro- 
vided for the different values of /nl in Table 1. 
We note that since the star formation rate densities are similar for 
/nl ^5 100, the best-fit parameters obtained do not differ signifi- 
cantly with respect to the Gaussian case. Also in the most extreme 
case /nl=1000 the GRB LF best-fit parameters are still consis- 
tent with those obtained in the Gaussian cosmology. However, in 
this case, the normalization k and the evolution parameter 8 are 
affected by the different shape of the cosmic SFR. This was in- 
deed expected: because of the enhanced SFR at high redshift in the 
/ NL =1000 model, a slightly lower evolution is required to repro- 
duce the observed redshift distribution of the complete sample of 
bright Swift GRBs and, consequently, also a different normaliza- 
tion is found. 

The LF of popIII GRBs is completely unknown. To compute their 
rate, we follow Campisi et al. (2011) and assume that popIII 
GRBs can be described by equation (7) with £ = 1.7 and L* = 
10 54 ergs -1 constant in redshift (i.e. L*,o = 10 54 ergs -1 and 
(5 = 0) (Toma et al. 2011). The normalization is then obtained by 
imposing that none of the ~ 500 GRBs detected by Swift so far 
were powered by popIII star explosions. We checked that our re- 
sults do not change significantly when varying £ between 1.5 and 2 
and log(L*,o/erg s -1 ) between 53 and 55. 

3.2 Validation from the Swift redshift distribution 

Before calculating the GRB rate expected for different cosmolo- 
gies, we test the validity of our theoretical model by means of the 
Swift data. We remind that, as of today, the Swift instrument has 
detected 604 GRBs in a lifespan of about 7 years, and the redshift 
complete (sub-)samples that have been extracted so far have only 
several tens of data points (see Perley et al. 2009; Greiner et al. 
201 1 ; Salvaterra et al. 2012; Hjorth et al. 2012). 
Fig. 2 reports the redshift evolution of all models expected at 
the Swift sensitivity, corresponding to a peak photon flux of 
0.4 ph s -1 cm -2 in the [15, 150] keV band. The Swift field of view 
of Ail s = 1-4 sr has been assumed. If we compare the /nl=0 and 
1000 cases, it is evident that significant differences arise at z > 6, 
where N changes by at least factor of ~ 2. At lower redshift the two 
distributions are very similar and possible differences fall within 
the uncertainties (shaded regions) on the evolution parameter, 8. In- 
deed, the upper and lower bounds of the shaded regions correspond 
to the evolution obtained by fitting the complete Swift sample with 
the maximum and minimum values of 8 as quoted in Table 1. We 
note that in principle an instrument like Swift can distinguish be- 
tween a Gaussian and a highly non-Gaussian model simply on the 
basis of the rate of GRB detections at high z. 
The four confirmed detections at z > 6 (GRB 050904 at z — 6.3 
by Kawai et al. 2006; GRB 080913 at z = 6.7 by Greiner et al. 
2009; GRB 090423 at z = 8.2 by Salvaterra et al. 2009,Tanvir et 
al. 2009; GRB 090429B at z ~ 9.4 by Cucchiara et al. 2011) cor- 
respond to a rate of JV(6) = 0.57 ± 0.28 GRBs per year, derived 
by using the entire timespan of Swift (~ 7 years). At face value, 
this is fully consistent with the predictions that our model provides 
for the Gaussian case. Moreover, since the GRB redshift distribu- 
tion for mildly non-Gaussian models does not differ significantly in 



© RAS, MNRAS 000, 1-12 



6 



U. Maio et al. 



Model 


log 10 (fc [M©" 1 ] ) 

Oil) V L W J / 


^*,0,51 


5 


<5 


C-stat 


f NL =o 
f NL =io 

f NL =50 

f NL =ioo 

f NL =1000 


7 7n +0.09 
'■' u -0.06 
7 7n +0.09 
— ' • ' u -0.06 
7 7 n +0.09 
' ■ ' u -0.06 
7 7n +0.09 
' ■ ' u -0.06 
7 7 r+0.10 
' ■ ' J -0.06 


o.iT±g;?g 
0.18±°;? 
o.i8±8;?g 

o-i9±8:S 

0.36i«;« 


2 n4 +0 ' 15 
2 04 +0 - 15 

2 04+° 15 
z - u4 -0.11 

2 10+ ' 20 


2 kq+0.63 

9 t-7+0.62 

9 t-7+0.62 
z -°'-0.57 
9 t-o+0.62 


29 
29 
29 
30 
30 



Table 1. Best-fit values and 1-cr errors for the free parameters of the GRB model, computed for the different cosmologies. From left to right the columns refer 
to the value of: /nl! me GRB normalization in [Mq -1 ], k; the characteristic luminosity at z = in [10 51 erg s — 1 ], L* o,5l! me slope parameter of the 
GRB LF, £; the redshift evolution parameter of the characteristic luminosity in the GRB LF, <5; the total C-stat value (i.e., the sum of the C-stat values obtained 
from the fit of the BATSE and Swift dataset) - for more details see Salvaterra et al. (2012). The total number of data points used to perform the fit is 33. 



the redshift range probed by Swift, the observed high-z rate is also 
consistent with any non-Gaussian model with a positive but smaller 
than - 100 /nl- 

However, we have to remind that the observed value for N(z) at 
z — 6 is a lower limit for the rate of high- z GRB detections with 
Swift, since some bursts at z > 6 could be hidden among the large 
fraction (~ 2/3) of GRBs for which the redshift has not been mea- 
sured. For this reason, the previous constraint seems to rule out 
very negative values of /nl that would fall below the aforemen- 
tioned limit. 

A strong upper limit of ^ 14% on the fraction of z > 6 GRBs 
detected by Swift has been recently determined by Jakobsson et al. 
(2012) 1 . Considering the 604 GRBs constituting the current Swift 
sample, this corresponds to at most 85 GRBs at z > 6, and to a rate 
of AT(6) < (12 ± 1) yr _1 (as this value is quite large, we do not 
show it in Fig. 2). 

The /nl=1000 case and larger values, already excluded by CMB 
analyses (Komatsu et al. 2011), lie off the observed rate by one or- 
der of magnitude or more. 

The previous considerations are based on the data point for the 
GRB rate at z ~ 6, but higher- z data, and larger, redshift complete 
samples (e.g. Perley et al. 2009; Greiner et al. 201 1; Salvaterra et al. 
2012; Hjorth et al. 2012) together with a better knowledge of the 
GRB luminosity function and of its redshift evolution, are needed 
in order to reduce error bars and to put tighter constraints on the 
amount of primordial non-Gaussianity on the basis of the observed 
GRB redshift distribution. More precisely, in order to discriminate, 
at redshift z — 6, between e.g. the /nl= and the /nl= 1000 
cases (whose number rates at z = 6 differ by a factor of ~ 2), with 
a 3cr confidence level, one should have a redshift complete sam- 
ple of roughly 800 GRBs (Poissonian errors have been assumed, 
and a current Swift number rate of N(0) = 88 yr _1 , according to 
Fig. 2). Such a large sample would also allow us to constrain the 
GRB LF quite accurately and then strongly reduce the error bars. 
A confidence level of la would require a smaller redshift complete 
sample of about 200 GRBs. 



3.3 Resulting GRB rates from the simulated cosmologies 

After the validation of the model (Sect. 3.2) performed by exploit- 
ing the calibration of the GRB rates against Swift data (Sect. 3.1.1 



1 We stress that the upper limit at z = 6 of 14% suggested by Jakobs- 
son et al. (2012) refers to a subsample of 69 Swift bursts and is obtained 
by assuming that the GRBs that could not be identified as low-redshift are 
actually at z > 6. Thus, the value of 14% must be taken as a very strong 
upper limit. 



and 3.1.2), we now apply it to different non-Gaussian cosmolo- 
gies to derive count predictions as a function of z. Note that we 
will assume the normalization derived from the simulations with 
100 Mpc/A-side boxes for the small 0.5 Mpc//?-side boxes, as well. 
In fact, the latter are not run down to z ~ and thus can not be 
used for normalization purposes. 

We stress that the following results are obtained by assuming an 
ideal instrument, that is able to detect all the GRBs at high red- 
shift. This is important to note, because, independently of the over- 
all normalization, the main effects of primordial non-Gaussianities 
on GRBs are expected to be originated by the differences shown in 
the redshift evolution of the SFRs (see Maio & Iannuzzi 201 1) and, 
hence, in the different GRB rates in the various models. 

3.3.1 Evolution of the GRB rates 

In Fig. 3 we plot the GRB rate, R, for the small 0.5 Mpc//i-size 
boxes (left panels) and the large 100 Mpc/A-size boxes (right pan- 
els). In the top panels, we show the redshift evolution for all the 
/nl scenarios considered, while in the bottom panels we focus on 
the relative contribution of the popIII GRB rate (Rm) to the to- 
tal rate, that is widely dominated by popII-I stellar generations, at 
z < 20. 

Besides small differences for the onset times of star formation, due 
to the different resolutions of the 0.5 and 100 Mpc//i-side boxes (see 
details on resolution issues in Maio et al. 2010; Maio & Iannuzzi 
201 1), in both small and large volumes the effects due to the pres- 
ence of primordial non-Gaussianities are visible at z > 10 — 15, 
while the rates eventually converge at later times, when feedback 
mechanisms start dominating the gas behaviour and the resulting 
star formation. 

In the small boxes, deviations from the Gaussian case are evident at 
earlier times, because these simulations can sample the very small 
primordial mini-haloes, which are extremely sensitive to the un- 
derlying matter distribution (top-left panel). As a consequence, star 
formation is resolved already at very high redshift and GRB rates 
of the order of ~ 10 -6 yr _1 sr _1 are expected at redshifts as high 
as z ~ 23 for /nl— 1000, and z ~ 19 - 20 for /nl— - 100. 
Similar values are reached in the large 100 Mpc/A-side boxes only 
at z ~ 20 for /nl— 1000, and z ~ 15 for /nl— - 100. 
These trends are valid for both popII-I and popIII regimes, even 
though the latter is usually negligible, predicting popIII GRB rates, 
Jim, that, following the behaviour of the popIII SFRs, drop by two 
orders of magnitude (bottom-left panel). 

The larger boxes miss the very small primordial haloes because 
of lack of resolution, but can sample much larger scales, showing 
that the effects of primordial non-Gaussianity can still be present 
at redshift z ~ 5 — 10 (top-right panel), i.e. for the whole first 
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Figure 3. Top panels: The expected popII-I and popIII GRB rates, R, in the 0.5 Mpc/ft-side boxes (left), and in the 100 Mpc/ft-side boxes (right), for 
models with different primordial non-Gaussianity: /nl= (solid black line), 10 (dotted cyan line), 50 (short-dashed green line), 100 (long-dashed red line), 
and 1000 (dot-dashed blue line). Bottom panels: The corresponding relative contributions of the popIII GRB rates, Rm, to the total rate, R, for the different 
cosmological models. 



billion years of the Universe, when the GRB rates should be only 
one or two orders of magnitude smaller than at present time. Also 
in these boxes, significant differences in the GRB rates are found 
only between the /nl = and /nl = 1000 scenarios. This holds 
for the corresponding popIII contributions (bottom-right panel), as 
well, and is consistent with what found in the smaller boxes, and 
with the converging behaviours at redshift below ~ 6. 

3.3.2 Comparison of the Gaussian and non-Gaussian models 

To directly compare and isolate non-Gaussian effects, in the up- 
per panels of Fig. 4, we plot the ratios between the results for the 
different non-Gaussian cases and the Gaussian model (/ NL =0), for 
both the 0.5 Mpc//i(left panels) and 100 Mpc/h(right panels) side 
boxes. Effects of large non-Gaussianities (/nl = 1000) are very 
well visible at almost any redshift with a rate that is boosted by 
about ~ 3 orders of magnitude in all boxes, at early epochs. This 
is due to the fact that in such models, over-densities are heavily 
biased to larger values, and, therefore, induce an earlier onset of 
star formation. More precisely, small scales (left panels) seem to 
depend very tightly on the underlying matter distribution, with en- 
hancements of the GRB rate at z ~ 20 by a factor of ~ 10 3 , 10, 3, 
and a few per cent, for /nl = 1000, 100, 50, and 10 respectively. 
In the latter three cases, feedback mechanisms from on-going star 
formation are able to reshuffle the gas and drive its hydrodynamical 
behaviour. As a consequence, the non-Gaussian effects below Mpc 
scales are almost washed out by redshift z ~ 15. The highly non- 
Gaussian case (/nl=1000), instead, shows more prolonged effects, 
with variations by a factor ~ 10 at z ~ 15, and a factor of ~ 2 at 
z Z 10. 

On larger scales (right panels), the ratios are similar for the various 
non-Gaussian cases, with corresponding delays toward lower red- 
shift in the 1ow-/nl scenarios. 



As a conclusion, we can state that the presence of primordial non- 
Gaussianities in the density fluctuations enhances early GRB rates 
and has effects up to z ~ 10 on Mpc scales, and at least z ~ 5 on 
much larger scales. 

To check whether different stellar populations can have different 
contributions, we can compare the corresponding ratios for popIII 
GRB rates only. The bottom panels in Fig. 4 readily demonstrate 
that the popIII GRB rates are less indicative of primordial non- 
Gaussianities, and less powerful in discriminating different /nl 
scenarios, mostly for /nl % 100. The fundamental reason is that 
the popIII contribution to the SFRs is very noisy due to the short 
life-times involved (Maio et al. 2010), and thus also the correspond- 
ing contribution to the GRB rates present more irregularities com- 
pared to the total (popII-I) GRB rates. 

Finally, we checked that uncertainties in the unknown popIII IMF 2 , 
in the Z cr u value, and in the stellar yields would not lead to sig- 
nificant differences for the previous results (see also related dis- 
cussions in Maio et al. 2010; Maio & Iannuzzi 2011). Similarly, 
changes in the popIII GRB efficiency do not alter the relative ef- 
fects of non-Gaussianities, since they would correspond just to a 
different overall normalization. 



4 DISCUSSION AND CONCLUSIONS 

In this work we have discussed the possibility of using GRBs as 
possible probe of the presence of primordial non-Gaussianities in 
the density field. This has been done using the outputs of two sets of 
N-body, hydrodynamic, chemistry simulations presented in Maio & 
Iannuzzi (201 1) (as also briefly described in Sect. 2). 

2 Here we adopted alternatively, as an extreme case, a Salpeter-like popIII 
IMF. 
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Figure 4. Top panels: The ratio between the popII-I GRB rates predicted for the different non-Gaussian models and the Gaussian case. Results for the 
0.5 Mpc/ ft-side boxes and the 100 Mpc/ ft-side boxes are shown in the left and right panels, respectively. Different lines refer to /nl= (solid black line), 
10 (dotted cyan line), 50 (short-dashed green line), 100 (long-dashed red line), 1000 (dot-dashed blue line). Bottom panels: The corresponding ratios for the 
popIII GRB rates. 



Besides gravity and hydrodynamics, the runs include radiative 
gas cooling both from molecules and atomic (resonant and fine- 
structure) transitions, star formation, UV background, wind feed- 
back, and chemistry evolution for various metal species, for both 
population III and population II-I stellar generations. 
Assuming that long 7-ray bursts are fair tracers of star formation 
(as suggested by e.g. Nuza et al. 2007; Lapi et al. 2008; Levesque 
et al. 2010; Campisi et al. 2011; Mannucci et al. 2011; Sanders 
et al. 2012; Michalowski et al. 2012), we propose to use them as 
probes of the underlying matter distribution at high redshift, when 
the possible presence of non-Gaussianity would have the strongest 
visible effects on the baryon evolution. 

By validating our calculations of the GRB rate against Swift data, 
we are able to exclude from the non-linearity parameter space very 
negative values for /nl (consistently with independent results from 
CMB data, e.g. Komatsu et al. 2011). 

When applying our model to different non-Gaussian scenarios, we 
find that already at 2 > 6 cosmologies with large /nl values 
present distinctive characteristics compared to those of the Gaus- 
sian case, independently from the errors on the LF parameter es- 
timates. Both on large and small scales, at very early times (2 ~ 
15 — 20) the boost in the rate due to non-Gaussianities is ~ 2 — 3 
orders of magnitudes for /nl=1000, and up to a factor of ~ 10 for 
/ NL =100. Differences of a factor of ~ 2 are still visible for milder 
values (/nl~ 50). However, while at small scales we find quick 
converging trends at lower redshift (2 ~ 9), determined by the lo- 
cally on-going star formation and feedback episodes, larger-scale 
volumes sample bigger objects and thus can retain memory of the 
primordial matter distribution even at 2 ~ 5 — 10. 
These effects are particularly evident on the total GRB rate, that is 
largely dominated by popII-I stars, while the result for the popIII 
GRB rate is noisier, mostly for /nl~ — 100, as a consequence 
of the corresponding, short-lived, popIII star forming regime (Maio 



et al. 2010). 

Additional changes in the popIII IMF, yields, Z cr it, or the overall 
normalization of the GRB rates will not alter these findings (see 
also more discussion in Maio et al. 2010). 

We have to recall that, when estimating the level of primordial 
non-Gaussianity, some difficulties come from the well known de- 
generacies of /nl with other factors, like cosmological parameters 
(e.g. the power spectrum normalization erg, or the equation-of-state 
parameter w), or with higher-order effects coming from baryonic 
matter evolution (e.g. supersonic bulk flows at early times; Tseli- 
akhovich & Hirata 2010; Maio et al. 201 1). 

We also warn the reader that the main assumption underlying our 
formalism is that GRBs are unbiased tracers of star formation 
(Fynbo et al. 2003, 2009; Stanek et al. 2006; Modjaz et al. 2008; 
Levesque et al. 2010; Levesque et al. 2010). Despite this has been 
recently supported by several works (see above), arguments for the 
existence of some possible biases exist in the literature, in particu- 
lar linked to metallicity selection of the host galaxies (e.g. Langer 
& Norman 2006). Such effects could alter the intrinsic redshift dis- 
tribution of GRBs (e.g. Natarajan et al. 2005; Langer & Norman 
2006; Salvaterra & Chincarini 2007; Cao et al. 2011; Salvaterra 
et al. 2012). However, we do not expect this to have significant 
impacts on the estimated trends of primordial non-Gaussianities. 
Indeed, any metallicity bias for the GRB formation is not supposed 
to be too strong, i.e. possible metallicity thresholds for the GRB 
progenitor stars can not be much lower than ~ 0.3 Zq (see Camp- 
isi etal. 2011). 

As shown in this paper, the differences in the GRB rate induced by 
non-Gaussianities are expected to be significant at very high red- 
shift. At 2 > 6 most of the galaxies (Salvaterra et al. 2011) and, 
in particular, most of the GRB progenitors (Salvaterra et al. 2012) 
have metallicities below this threshold (see also detailed studies 
in Maio et al. 2010). These studies find that only a small fraction 
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(< 5%) of galaxies at z = 6 has got a metallicity Z ^ 0.3 Z©, 
rapidly decreasing at higher redshift. Therefore, at least at these 
early times, our assumption of GRBs as fair tracers of the cosmic 
SFR is quite solid. Furthermore, we checked that the difference 
among Gaussian and non-Gausssian models remains unchanged 
when galaxies with metallicities larger than 0.3 Zq were excluded 
from our analyses. 

We note that estimates of non-Gaussianities (e.g. Komatsu et al. 
2011) based on cosmic microwave background and large-scale- 
structure data seem to support positive /nl values up to ~ 100. 
This implies that at early epochs we expect an enhancement of the 
GRB rate up to a factor of 10 with respect to the standard Gaussian 
case. 

We stress that the existence of GRBs at such high redshift is not 
unlikely, as they are tightly linked to star formation episodes. In 
principle, they could be observable thanks to their large intrin- 
sic luminosity and longer time dilution of the afterglow. None the 
less, from an observational point of view, detections of GRB after- 
glows at very high redshift are complicated by Lyman-a absorption 
from inter-galactic gas. In fact, for bursts at z ^ 15, as the ones 
we are interested in here, no flux can be detected in photometric 
bands bluer than the K band (at ~ 2.2 fim). At z > 18, where 
the largest differences between Gaussian and mildly non-Gaussian 
models are expected, observations in the infra-red band are needed. 
Since follow-up observations of GRB afterglow are generally car- 
ried out in optical-NIR bands, extreme high-z GRBs can be missed. 
However, a small population of extremely dark GRBs (e.g. Greiner 
et al. 201 1), i.e. bursts for which the afterglow remains undetected 
in spite of early and deep K band observations, has been recently 
identified (D'Elia & Stratta 201 1). While the nature of these GRBs 
is still matter of debate and alternative explanations for their dark- 
ness do exist 3 , it is possible that these bursts (or at least one of 
them) are at z ^ 18. If confirmed, this could provide evidence 
in favor of a mildly positive non-Gaussian parameter (/nl in the 
range 10 — 100, see Fig. 3). Future detections of extremely dark 
GRBs (as the ones by D'Elia & Stratta 2011) at redshift z > 20 
and with a substantial rate, of at least ~ 10~ 6 yr -1 sr _1 , might 
be an indication of even bigger values for /nl. Naively speak- 
ing, a determination of the rate for such GRB would lead to about 
~ (0.1±0.1) yr _1 sr _1 , but one should also consider that the prob- 
ability of observing such event is almost as small as ~ 10 -3 , since 
this is a unique case out of the 604 Swift GRBs. In principle, this 
would imply positive /nl values, but with huge error bars. 
However, in order to draw more definitive conclusions and give 
more stringent constraints much larger high-z GRB complete sam- 
ples, currently not available in the literature, are required. 
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